High-accuracy calculation of black-body radiation shift in ^^'^Cs primary frequency 
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Black-body radiation (BBR) shift is an important systematic correction for the atomic frequency 
standards realizing the SI unit of time. Presently, there is a controversy over the value of the BBR 
shift for the primary ^'^^Cs standard. At room temperatures the values from various groups differ at 
3 X 10~^^ level, while the modern clocks are aiming at 10"^'' accuracies. We carry out high-precision 
relativistic many-body calculations of the BBR shift. For the BBR coefficient f3 a,t T — 3007^ we 
obtain /3 — — (1.708±0.006) x 10"^'', implying 6x 10~^^ fractional uncertainty. While in accord with 
the most accurate measurement, our 0.35%-accurate value is in a substantial, 10%, disagreement 
with recent semi-empirical calculations. We identify an oversight in those calculations. 

PACS numbers: 06.30.Ft, 32.10.Dk, 31.25.-v 
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Since 1967, the SI unit of time, the second, is defined 
as a duration of a certain number of periods of radiation 
corresponding to the transition between two hyperfine 
levels {F — A and F = 3) of the ground state of the 
^^■^Cs atom. In 1997, this definition has been amended, 
to specify that the above statement refers to the atom 
at a temperature of K (and at rest) 0. The atomic 
clocks are usually operated at room temperatures and 
the specific reference to T = K implies that a proper 
correction for the action of the thermal bath of photons 
on the atomic energy levels is explicitly introduced. For 
^^•^Cs clocks, it is an important systematic correction 0, 
as the resulting fractional correction to atomic frequency 
due to black-body radiation (BBR) at T = 300 is in the 
order of 10~^^. Moreover, presently, there is a contro- 
versy over the value of the BBR shift for the primary 
^^•^Cs standard. At T = 300 K the values from various 
groups 3, 4, 5, 6, 7] differ at 3x 10~^^ level, while modern 
Cs clocks aim at 10^^® accuracies @- 

The persistent discrepancies in the BBR shift have 
prompted the efforts at the US National Institute for 
Standards and Technology (NIST) on removing the sensi- 
tivity to BBR by operating the primary standard at cryo- 
genic temperatures 9j (BBR shift scales as T*). How- 
ever, because of the weight limitations, this direct ap- 
proach would be hardy feasible if next-generation atomic 
clock were to be operated at the International Space 
Station This ongoing controversy and implications 
for atomic time-keeping serve as motivations for our pa- 
per. Here we compute the ^^'^Cs BBR shift using high- 
accuracy relativistic many-body techniques of atomic 
structure. Our evaluated error in the BBR shift implies 
6 X 10^^^ fractional uncertainty in the clock frequency 
with the value of the BBR shift consistent with the most 
accurate (0.2%-accurate) measurement However, our 
0.35%-accurate value is in a substantial 10% disagree- 
ment with recent semi-empirical calculations 0, • We 
show that this discrepancy is due to contributions of the 
intermediate continuum states omitted in those calcula- 
tions. 

First, let us review underlying theory of the Cs BBR 



shift. BBR causes a weak oscillating perturbation of 
atomic energy levels. Conventionally, the leading term 
in the BBR contribution is parameterized as the frac- 
tional correction to the unperturbed clock frequency, 
1^0 = 9192 631770 Hz, 

where Tq = 300 K. Evaluation of the coefficient P is 
the goal of this work. This coefficient can be related 
to the scalar differential polarizability for the hyperfine 
manifold of the 6si/2 Cs ground state. Indeed, the char- 
acteristic thermal photon energy at room temperatures 
is much smaller then the atomic energies, so that the 
perturbation can be well described in the static limit. 
Moreover, contributions of electro-magnetic BBR multi- 
poles beyond electric dipoles, as well as retardation cor- 
rections, are highly suppressed 0|. Then the BBR shift 
of the energy level is given by (atomic units are used 
throughout, a is the fine-structure constant) 



— (a^)3rV(0), 



(1) 



where aF{0) is the static scalar electric-dipole polariz- 
ability of the hyperfine level F. The vector and tensor 
parts of the polarizability average out due to the isotropic 
nature of the BBR. 

The relation of the BBR shift to polarizability has 
been exploited in the most accurate measurement of 
the differential Stark shift However, recent di- 

rect temperature-dependent measurement of the BBR 
shift turned out to be differing by about two standard 
deviations from the indirect measurement. Namely this 
difference has stimulated the recent interest in the Cs 
BBR shift. 

The overall BBR shift of the clock frequency is the 
difference of the individual shifts for the two hyperfine 
states {F — A and F — 3) involved in the transition. 
While taking the difference, the traditional lowest-order 
polarizability of the Qsi/2 level cancels out and one needs 
to evaluate the third-order .F-dependent polarizability 
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ap {0). This contribution involves two El interactions 
y = — D • E with the external electric field and one hy- 
perfine 7?hfs coupling We parameterize the general 
correction to the energy as a sum of four diagrams, 



can be written more explicitly in terms of the F- 
independent diagrams as 



SE 



(3) _ 



= {F\Hi,f,RVRV\F) + {F\VRHmsRV\F) + 

{F\VRVRHms\F) - {F\Hms\F){F\VRRV\F) . 



where R = {H — Ey)~-^ is the resolvent operator, with 
H and Ey being the atomic Hamiltonian and the ground 
state energy, respectively. The four contributions will 
give rise to top, center, bottom, and residual contribu- 
tions discussed below. (Naming convention reflects rel- 
ative position of the hyperfine operator in the string of 
the three interactions). The top and bottom diagrams 
are equal due to hermicity considerations. Due to angu- 
lar selection rules only magnetic dipolar interaction re- 
mains; we write H^s — M ' where ^ is the nuclear 
magnetic moment and T*^^) = (^a ■ C'j°''(f)^ /(cr^) 

is the relevant relativistic coupling tensor [T^ . 

After angular and many-body reduction, the scalar po- 
larizability may be expressed as 

a^piO) = iv/(2/)(2/+l)(2/ + 2)|^" f 

(-1)™'" {2T + C + R) 

where gi is the nuclear gyromagnetic ratio, fi„ is the 
nuclear magneton, / = 7/2 is the nuclear spin, and jy = 
1/2 is the total angular momentum of the ground state. 
The F-independent sums are {\v) = |6si/2)) 

^ {v \\D\\ to) (m ||L>|| n) {n \\T^^^ \\ v) 



2jy + 1 



{Em — Ey){En — Ey) 
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The summation indexes to and n range over valence 
bound and continuum many-body states and also over 
single-particle core orbitals. With this convention, the 
above expressions subsume contributions from intermedi- 
ate valence and core-excited states and they also take into 
account so-called core- valence counter-terms |0. Selec- 
tion rules impose the following angular symmetries on the 
intermediate states: si/2 for \n), Pi/2,3/2 for |to) in the 
top diagram, Pi/2.3/2 for both \m) and \n) in the center 
diagram, and, finally, J31/2.3/2 for the \m) in the residual 
term. 

We will tabulate our results in terms of the 
conventional scalar Stark shift coefficient kg — 

-l/2(aP^^{0) - Jpl^{0)y For ^^^Cs this coefficient 
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giUn {2T + C + R) , 



(2) 



with the BBR coefficient (3 = -4/15 {an)^ T^/vq x kg. 

Numerical evaluation of the diagrams T,C, and R can 
be carried out either using the Dalgarno-Lewis method or 
by directly summing over individual intermediate states. 
Here we use the direct summation approach. This treat- 
ment is similar to the one used in high-accuracy cal- 
culations of atomic parity violation in ^'^^Cs HJ]. The 
main advantage of this method is that one could explic- 
itly exploit high-accuracy experimental data for energies, 
dipole-matrix elements, and hyperfine constants. When 
the accurate values are unknown, we use ah initio data 
of proper accuracy. In addition, this approach facilitates 
comparison with recent calculations |a, 0] , which also use 
the direct summation approach. 

The central technical issue arising in direct summation 
over a complete set of states is representation of the in- 
numerable spectrum of atomic states. For example, even 
without the continuum, the bound spectrum contains an 
infinite number of states. A powerful numerical method 
for reducing the infinite summations/integrations to a fi- 
nite number of contributions is the basis set technique. 
In particular, we employ the B-spline technique 15]. In 
this approach an atom is placed in a large spherical cav- 
ity and the single-particle Dirac-Hartree-Fock (DHF) or- 
bitals are expanded in terms of a finite set of B-splincs. 
The expansion coefficients are obtained by invoking vari- 
ational Galerkin principle. The resulting set of the single- 
'fjarticle orbitals is numerically complete and finite. The 
technique has a high numerical accuracy and we refer 
the reader to a review [l^ on numerous applications of 
B-splines for details. 

We use B-spline set with 70 splines of order 7 for each 
partial wave and constrain the orbitals to a cavity of 
radius i?cav = 220 a.u. This particular choice of i?cav 
ensures that the lowest-energy atomic orbitals are not 
perturbed by the cavity. In particular, all core and va- 
lence DHF orbitals with radial quantum numbers 1 — 12 
from the basis set produce energies and matrix elements 
in a close numerical agreement with the data from tra- 
ditional finite-difference DHF code. These low-energy 
orbitals will produce true many-body states for a cavity- 
unconstrained atom. 

To understand the relative role of various contribu- 
tions, we start by computing the Stark shift at the DHF 
level. We obtain: 2Tff^ _ n^is 



= 0.418, = 0.003, 



and -rmrr 



0.518, resulting in the Stark coefficient of 

fcDHF = -2.799 X IQ-io Hz/(V/m)2 [l3. It is clear that 
the top and residual terms dominate over the center di- 
agram. The bulk (99.8%) of the value of the residual 
term is accumulated due to the principal 6s — 6^1/2.3/2 
transitions. For the top term the saturation of the sum is 



3 



not as rapid, but still the dominant contributions come 
from the lowest-energy excitations: limiting the summa- 
tions to the first four excited states recovers only 68% 
of the total value for the top diagram. In addition, we 
find that core-excited states contribute only 0.04% to the 
final value. 

The above observations determine our strategy for 
more accurate calculations. We group the entire set 
of atomic states into the "main" low-lying-energy states 
(principal quantum numbers n <= 12) and remaining 
"tail" states. We will account for the contribution from 
the "main" states using high-accuracy experimental and 
ab initio values. The contribution from the "tail" will be 
obtained using either DHF or mixed approach. 

First, we describe the high-accuracy data used in 
our calculations. We need dipole and hyperfine ma- 
trix elements and energies. Experimental values for 
the dipole matrix elements for the following six transi- 
tions were taken from the literature (see compilations in 

Refs. Hill) 6Si/2 - 6pi/2,3/2,7Si/2 - 6^1/2,3/2, 7si/2 - 

7pi/2,3/2- Crucial to the accuracy of the present analy- 
sis were the matrix elements for the principal 6si/2 — 



6p 



1/2,3/2 



transitions. We have used 0.005%-accurate 



value for (6s 1/2 H^H 6^3/2) from Ref. [19j. The value 
for (6S1/2II-DII6P1/2) was obtained by using the above 
6si/2 — 6p3/2 matrix element and 0.03%-accurate mea- 
sured ratio (22j of these matrix elements. These six ex- 
perimental matrix elements were supplemented by 92 
values {nsi/2 — n'pi/2.3/2 values for n,n' — 6 — 12 ) 
from high-accuracy ab initio calculations. We employ 
the relativistic linearized coupled-cluster singles-doubles 
(LCCSD) method. The underlying formalism, imple- 
mentation, and results for alkali atoms are described in 
Ref. ^3 . For dipole matrix elements the accuracy of the 
ab initio LCCSD method is a few 0.1%. 

As to the high-accuracy values of the matrix el- 
ements of the hyperfine coupling, the diagonal ma- 
trix elements of the T'-^-' tensor are directly re- 
lated to the conventional hyperfine constants: A = 
gifijj, [(2j.)/(2j„ + l)/(2j, + 2)f^ {v\\T(^mv). Wc 
have used the compilation of hyperfine constants from 
Ref. [231 for the "main" n = 6 — 12 states. Off-diagonal 
matrix elements between the s-states were evaluated us- 
ing the geometric-mean formula 



nsi/2 



Si/2 



{|< 



'^Sl/2 



nsi/2)(n'si/2 



^'si/2> } 



1/2 



This expression has been shown to hold to about 10^"^ in 
Ref. [23 (notice that the radiative corrections would start 
playing a role at a few 0.1% as well). If we had 6 < n < 12 
in the above expression, then the experimental value is 
used for its corresponding diagonal element on the right. 
If we also had 6 < n' < 12, then that experimental value 
is also used for the corresponding diagonal element on 
the right, otherwise the DHF value is taken, (n and n' 
can be obviously interchanged in this prescription.) This 



mixed approach has allowed us to uniformly improve the 
accuracy of the calculations. Indeed, in the numerically 
important top term, the hyperfine matrix elements come 
in the combination (nsx/2 ||^*'^^|| 651/2)- As n grows, the 
correlations become less important, so the dominant cor- 
relation correction comes from the 6si/2 state. Using 
the described mixed approach allows us to account for 
this dominant correlation. The geometric-mean formula 
holds only for the s states. For the off-diagonal matrix 
elements between various combinations of 6^1/2,3/2 ^^^id 
^Pi/2,3/2 (" = 6 — 9) states we employed a modifica- 
tion of the LCCSD method augmented by perturbative 
treatment of the valence triple excitations (LCCSDpvT 
method) The accuracy of these matrix elements is 

a few %. As these matrix elements enter the relatively 
small center term, the effect on the overall theoretical 
error is negligible. 

Finally, we used experimental energy values from the 
NIST tabulation ^23], for states with principle quantum 
number rt = 6 — 12, and the DHF values otherwise. 

With the described set, we report our final result for 
the scalar Stark coefficient and relative blackbody radia- 
tion shift at 300K to be 



-(2.268 ± 0.008) X 10"^° Hz/(V/m)2 , 
-(1.708 ±0.006) X 10"^^ 



(3) 



The values of the individual diagrams are |^ = 0.449, 



^ = -0.002, and ^ = 0.553 

ks ' fc - 



ks 



When comparing with 
the DHF values, the most substantial modification due 
to correlations is in the center term, which changes the 
sign. Fortunately, this term is relatively small, and this 
extreme change does not substantially affect the final re- 
sult. 

The overall uncertainty of these results was determined 
from the uncertainties of the individual matrix elements 
and energy values used in their computation. Standard 
uncertainty analysis was done throughout all mathemat- 
ical operations. For energy values taken from NIST, the 
uncertainty is assumed negligible. For all other experi- 
mental values, the reported uncertainty is used. For ab 
initio matrix elements (DHF, LCCSD, or LCCSDpvT) 
we assigned an assumed uncertainty. These assumed un- 
certainties were based on comparison between calculated 
and high-accuracy experimental values. This resulted in 
a relative uncertainty for both the scalar Stark coefficient 
and the BBR shift of 0.35%. We have performed several 
consistency checks, e.g., replacing experimental matrix 
elements and energies by ab initio LCCSD values or by 
replacing the DHF values for states with n = 13 — 27 with 
the LCCSD values. The final result was stable to such 
modifications within the stated uncertainty in Eq.(|3J)- 
These tests provide us with additional confidence with 
respect to our standard error analysis based on errors of 
used experimental values. It is also worth noting that the 
present calculation does not include radiative corrections 
which may contribute at the level of a few 0.1% (some 
radiative corrections, e.g., vacuum polarization, are ab- 
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TABLE I: Values of ks in 10-^°Hz/(V/m)^ 



References 

theory -1.97± 0.09 Ref. [6] 

theory -2.06± 0.01 Ref. [T] 

expt. -2.05± 0.04 Ref. \£ 

expt. -2.271± 0.004 Ref. [3j 

theory -2.268± 0.008 present 



sorbed in our final value already as we use experimental 
hyperfine constants). 

A comparison with recent theoretical and experimen- 
tal work is presented in Table |I| While agreein_g with 
the most accurate measurement by Simon et al. !3], our 
results are in substantial disagreement with the recent 
calculations |^ . The principal differences between the 
present work and these calculations are: (i) more sophis- 
ticated treatment of correlations, and (ii) rigorous sum- 
mation over the complete set of intermediate states in 
perturbative expressions. As discussed above, we used 
the numerically complete basis-set approach which ap- 
proximates Rydberg states and continuum with quasi- 
spectrum. By contrast, in Ref. ,6|, the summations were 
truncated at n = 9, and in Ref. |2] at n = 18; neither 
work includes continuum. To illuminate the importance 
of the omitted contributions we truncate our summations 



at n = 12. The resulting value deviates from our final ks 
result by 0.29 x 10^^° Hz/(V/m)2. This large 10% "con- 
tinuum correction" brings the values from Refs. 0, 0] 
into essential agreement with our result. The fact that 
continuum needs to be included is hardly surprising, as, 
for example, about 20% of the textbook polarizability of 
the ground state of the hydrogen atom comes from the 
continuum states. 

To conclude, here we have reported results of rela- 
tivistic many-body calculations of the BBR shift, one 
of the leading systematic correction in ^^'^Cs frequency 
standard and a subject of the recent controversy. Our 
0.35%-accurate result re-validates high-precision Stark 
shift measurements 0] . Our work also clarifies the origin 
of the reported discrepancy between that measurement 
and recent calculations |^ . 

We would hke to thank H. Gould, M.-A. Bouchiat, 
and W. Itano for discussions. Work of K. B and A.D 
was supported in part by NIST Precision Measurement 
Grant, and National Science Foundation. Work of U.I.S 
was supported in part by DOE-NNSA/NV Gooperative 
Agreement DE-FC52-01NV14050. 

Note. While completing writing this manuscript, we 
have learned of another accurate many-body calculation 
of the BBR shift in '^^^Gs clock H^. Their result, = 
2.26 X 10^^^ ± l%Hz/(V/m)^, is in agreement with our 
more accurate value. 
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